Molecular dynamics simulations of the evaporation of particle- laden droplets 
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We use molecular dynamics simulations to study the evaporation of particle-laden droplets on a 
heated surface. The droplets are composed of a Lennard- Jones fluid containing rigid particles which 
are spherical sections of an atomic lattice, and heating is controlled through the temperature of an 
atomistic substrate. We observe that sufficiently large (but still nano-sized) particle-laden drops 
exhibit contact line pinning, measure the outward fluid flow field which advects particle to the drop 
rim, and find that the structure of the resulting aggregate varies with inter-particle interactions. In 
addition, the profile of the evaporative fluid flux is measured with and without particles present, 
and is also found to be in qualitative agreement with earlier theory. The compatibility of simple 
nanoscale calculations and micron-scale experiments indicates that molecular simulation may be 
used to predict aggregate structure in evaporative growth processes. 
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The evaporation of a sessile droplet on a hot surface is 
a key problem in fluid mechanics, relevant both to theo- 
retical issues in heat transfer and to practical questions 
in materials processing. The evaporation of a particle- 
laden droplet raises the additional issue of the structure 
of the resulting solid aggregate, and, going further, of- 
fers the possibility of controlling this structure by means 
of anisotropic (e.^., Janus) surface properties [l^. A fa- 
miliar and paradigmatic example of this process occurs 
in coffee stains, where the residue of evaporated droplets 
takes the form of a ring-like deposit of grains at the rim. 
Experiments by Deegan and collaborators [2| focused at- 
tention on this "coffee ring problem" several years ago, 
and subsequent work [3|-l6[ estabhshed the ubiquity of 
the process, while numerous theoretical studies have ad- 
dressed the dynamics |7H11|. A complete understand- 
ing of the problem is not yet available however: ex- 
periments cannot measure everything in a small, time- 
dependent, multiphase droplet, while most theoretical 
treatments require approximations to deal with an evap- 
orating particle-laden drop. 

In this paper we use molecular dynamics (MD) simula- 
tions to simulate the evaporation of droplets containing 
colloidal particles, having either uniform or Janus-like 
surface properties. One goal is to test whether the phe- 
nomena found in micron-sized particle systems persist 
down to nanometer scales; in this way we hope to extend 
the size range in which controlled aggregate structures 
may be produced by droplet evaporation. A second goal 
is to test the validity of some of the underpinnings of the 
theoretical analyses used in the problem. Since MD simu- 
lations provide detailed atomic-scale information, we can 
measure concentration, temperature and fluid flow, even 
during the rapid heterogeneous processes occurring in 
evaporation. The difficulties of applying uncertain con- 
stitutive relations are absent, although replaced to some 
degree by the problem of extracting a robust signal from 
a relatively small sample in a fluctuating environment. 
More generally, our goal is to establish the ability of these 



relatively basic simulations of moderate scale systems to 
predict phenomena occurring in droplet evaporation and 
guide experimental investigations. 

The simulations use standard molecular dynamics 
(MD) techniques 12l4l4j and generic interactions of 
Lennard- Jones form, 



V{n,)=Ae 
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For simplicity we assume that all fluid and particle atoms 
have the same interaction potential, along with the same 
mass m and approximate diameter <j. The calculations 
are nondimensionalized using e, <j and m as energy, 
length and mass scales, respectively, and the resulting 
time scale is r = cF{m/eY^'^ . Typical numerical values 
are a ~ 0.3nm, r ~ 2ps and e ~ 120A:b, where ks is 
Boltzmann's constant, and temperatures are measured 
in units of e/kB- The fluid atoms in liquid or vapor obey 
ordinary Newtonian dynamics with the force arising from 
the interaction with other atoms (within a cutoff radius 
of 2.5cr), using a predictor-corrector method for the in- 
tegration of the equations of motion. The particles are 
spherical sections of an atomic fee lattice containing all 
atoms within a certain radius of a center; here the atomic 
density is 0.8, the radius is 2 and the particles contain 
32 atoms. These are Janus particles, implemented by 
making the atoms in only one hemisphere of a particle at- 
tractive to other particle atoms; all particle atoms attract 
the fluid and wall in the same way. The particles move 
as rigid bodies, where the net force and torque on each 
particle is computed by summing the interatomic forces 
between its atoms and the neighboring fluid atoms, and 
the motion is given by the Newton and Euler equations. 
Quaternion variables are used to describe the particle 
orientations [12|. The solid substrate is made of atoms 
coupled to fee lattice sites by a linear spring of stiffness 
lOOe/cr^. 

Initially the drop consists of a hemispherical cap of 
72,236 fluid atoms placed above a (monolayer) solid wall. 
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puted by dividing the simulation domain into concentric 
circular rings centered along the vertical axis of the drop, 
and simply averaging the atomic velocities in each ring. 
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FIG. 1. (Color online) Evaporation of a particle- laden droplet: 
fluid and solid atoms are shown as cyan and red dots respec- 
tively and the Janus particles are circles whose two sides are 
red and blue. Top to bottom: times 50, 500 and lOOOr. 



FIG. 2. Velocity field in cylindrical coordinates in an evapo- 
rating particle-laden droplet at early and late times; the solid 
occupies the region z < 0. 

The origin of the flow that drives particles to the rim 
is, it is believed [2|, contact line pinning coupled to the 
fact that the evaporative flux is largest at the edges of 
the drop. The liquid must supply this flux as the droplet 
shrinks down, and the geometry of the situation requires 
a strong outward flow field, as seen in Fig. [2l More pre- 
cisely, for drops of pure liquid evaporating with a fixed 



contact line, one can show 
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Here, 6 is the contact angle between the droplet and the 
wall, which varies as the droplet evaporates. The key 
point is that for < 7r/2 the (mathematical) vapor flux 
diverges at the edge of the droplet. Of course, there is no 
real singularity in a physical problem, and one expects j 
to be cut off at a small (molecular) scale. 



119 particles of 32 atoms each are centered at random 
positions within the cap, and the entire system is equili- 
brated by gradually raising the temperature from 0.5 to 
1.2 using velocity rescaling, following which the substrate 
is maintained at this temperature while the fluid tem- 
perature is allowed to vary. The liquid expands slightly 
during the temperature ramp and then, as seen in Fig. (U 
the drop shrinks monotonically as it emits vapor, and 
eventually disappears due to evaporation. Evaporating 
fluid atoms which leave the simulation box are simply 
deleted. Meanwhile the particles settle towards the sub- 
strate and subsequently are advected to the rim of the 
droplet and deposit there. The contact line itself remains 
pinned. The connection between liquid and particle mo- 
tion is indicated by the velocity field shown in Fig. O the 
fluid moves downward over most of the drop and radially 
outward near the substrate. This velocity field is com- 
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FIG. 3. (Color online) Evaporative vapor flux vs. scaled 
radius for particle laden-drops and, in the inset, for pure fluid 
drops. The three curves refer to times 300r (x, green), 400t 
(A, blue) and 500t (D, orange). 



We have measured the evaporative flux both for pure 
hquid and particle-laden drops, and obtained Fig. [3l 
which shows j(r, t) at the three successive times indi- 
cated. Each plot is an average over a short (lOr) interval 
centered at 300, 400 and 500r. Averaging is necessary 
to smooth the fluctuations, and longer averaging periods 
would be more effective, but the drop shape changes too 
much over longer intervals. The entire evaporation pro- 
cess lasts for about lOOOr (roughly 2 ns) for the Janus 
case, and slightly longer for pure fluid evaporation, but 
the data at later times involves fewer and fewer evapo- 
rating atoms and is too noisy for analysis. The ordinate 
in Fig. [3] is radial position divided by the current drop ra- 
dius, and varies between and 1. While the radius is con- 
stant (except for fluctuations) in the particle-laden case, 
the radius of the pure fluid drop decreases with time. We 
see that indeed the vapor current tends to diverge at the 
rim and the divergence increases with time. This behav- 
ior is qualitatively consistent with Eq. [21 since a stronger 
divergence requires a smaller contact angle, and the sim- 
ulations show this trend. However, the formula fails to 
quantitatively describe the simulations. The exponents 
A(^) obtained in fitting the pure fluid data are 0.48, 0.55 
and 0.60 at times 300, 400 and 500r, respectively, but 
Eq. [2]does not permit values A > 1/2. Some possible ex- 
planations for this discrepancy are deviations from bulk 
continuum behavior in nanosized regions, difficulties in 
measuring the flux accurately in directions nearly paral- 
lel to the surface, and inadequate ensemble averaging in 
obtaining the data. The flux is slightly more singular for 
the particle-laden drop: the fitted values of A are 0.53, 
0.58 and 0.62 for the same three times. 

Not all evaporating nanodrops behave in this way. One 
important requirement for obtaining a deposit at the rim 
is that the drop be large enough: if the drop is too small 
it evaporates, or at least decays into a thin pancake, be- 
fore the flow field is established and the particles are able 
to move to the rim. This behavior was first observed ex- 
perimentally by Shen et al. [17| , and we have reproduced 
it in simulations. The drop shown in Fig. [T] has a radius 
of about 20 nm, whereas in similar simulations for drops 
whose initial radius is 5 or 10 nm we see that the par- 
ticles deposit roughly uniformly over the drops interior; 
see Fig. HI 

A second constraint needed for deposition to occur 
at the drop rim is that the liquid must have adequate 
thermal contact with the solid to set up the flow field 
seen above [18[ . We have investigated this issue by vary- 
ing the interaction between the liquid and the wall: in 
the Lennard Jones potential (Eq. [1]) acting between fluid 
and wall atoms, we varied the coefficient of the attrac- 
tive r~^ term between (pure short-distance repulsion - 
hydrophobic wall) and 1 (standard strength attraction - 
completely wetting wall) and observed the resulting solid 
pattern while measuring the thermal conductivity. The 
latter simulation involved a slab of liquid completely fill- 
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FIG. 4. (Color online) Particle deposit for evaporated drops 
of initial radius (left) 5 and (right) 10 nm. The color scheme 
is given in the caption to Fig. [1] 



ing the gap between two atomic walls held at different 
temperatures. As the attractive strength decreased from 
1 to the thermal conductivity decreased approximately 
linearly, ultimately by a factor of nearly 10. Correspond- 
ingly, as the wall attractive strength decreased the re- 
sulting pattern of solid particles varied from the rim de- 
posit shown above to a crystallized droplet. The effects 
of substrate thermal resistance and conductivity on evap- 
oration have been studied more systematically by Dunn 
et al 
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The detailed structure of the deposit is an impor- 
tant consideration in potential applications to evapora- 
tive self-assembly. We saw above that particles with 
Janus surface properties would, under the right con- 
ditions, form a rim deposit with chain-like structures. 
This behavior is confirmed by experimental observations 
of gold-capped sulfated-polystyrene Janus particle-laden 
droplets during drying [21]. The chaining behavior is a 
result of the attractive interaction between the gold caps 
of and the sulfated polystyrene half of the Janus particles. 

Additional experiments investigated the effect of par- 
ticle surface charge on the deposition behavior [2l|, |22 |. 
The presence of sulfate groups on the surface of the par- 
ticles results in an overall negative charge, which can be 
screened by the addition of salt. Sulfated polystyrene 
particles in distilled water form hexagonal close-packed, 
highly ordered layers at the rim, whereas formation of 
randomly packed particle layers is observed in 10 mM 
aqueous NaCl solution. These results motivated us to 
study the effect of charge on the deposit structure of uni- 
form particles. In Fig.[5l we show the deposits that result 
when the charge on a symmetric particle increases from 
to 4 and then to 8. In the simulations, a charge of the 
appropriate magnitude is placed at the center of the par- 
ticle, and a Coulomb interaction is applied. The neutral 
case resembles the random packing observed in the 10 
mM case. The charge 8 situation models the case where 
the particles carry a charge, i.e., deposition in deionized 
water. It is apparent from the bottom panel of Fig.[5]that 
the particles tend toward assembling at the rim and a 
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tration and stress fields within tiie droplet. In addition, 
we were able to show the existence of a minimum drop 
size for rim deposition, and verify the importance of ad- 
equate thermal coupling between liquid and solid. Fur- 
thermore, the simulations give the connection between 
particle interactions and deposit structure. The signifi- 
cance of these results is that simple simulations provide a 
viable method for both testing the theoretical underpin- 
nings of the process and for predicting the nature of the 
outcome - the structure of the resulting particle deposit. 
In this paper we have focused on the most important 
aspect of the continuum flow, the velocity field and the 
evaporative flux, but any other quantity which can be 
determined from atomic variables is equally accessible. 

This work was supported in part by NSF-CBET (CA- 
REER) 0644789. 
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FIG. 5. Effects of adding a charge to a plain (non- Janus) 
particle. Top to bottom: charge 0, 4 and 8. 



regular packing [235 ^ t)ut that the drop volume is not suf- 
ficient to enable long enough evaporation times to achieve 
close packing at the rim. 

We have shown that straightforward, medium-scale 
MD simulations can easily capture most of the salient 
features in the evaporation of particle-laden droplets. 
Aside from demonstrating that nano-scale and micron- 
scale systems behave in a similar way with regard to the 
behavior of the particles, we were able to measure con- 
tinuum fields such as velocity within the droplet, along 
with the profile of the evaporative flux, which drives the 
process. Using standard methods, we have also measured 
(but not reported here) the density, temperature, concen- 
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